Epidemic fronts in complex networks with metapopulation structure 
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Infection dynamics have been studied extensively on complex networks, yielding insight into the ef- 
fects of heterogeneity in contact patterns on disease spread. Somewhat separately, metapopulations 
have provided a paradigm for modeling systems with spatially extended and "patchy" organiza- 
tion. In this paper we expand on the use of multitype networks for combining these paradigms, 
such that simple contagion models can include complexity in the agent interactions and multi- 
scale structure. We first present a generalization of the Volz-Miller mean- field approximation for 
Susceptible-Infected-Recovered (SIR) dynamics on multitype networks. We then use this technique 
to study the special case of epidemic fronts propagating on a one-dimensional lattice of intercon- 
nected networks - representing a simple chain of coupled population centers - as a necessary first 
step in understanding how macro-scale disease spread depends on micro-scale topology. Using the 
formalism of front propagation into unstable states, we derive the effective transport coefficients of 
the linear spreading: asymptotic speed, characteristic wavelength, and diffusion coefficient for the 
leading edge of the pulled fronts, and analyze their dependence on the underlying graph structure. 
We also derive the epidemic threshold for the system and study the front profile for various network 
configurations. To our knowledge, this is the first such application of front propagation concepts to 
random network models. 

PACS numbers: 64.60.aq, 87.19.X-, 87.23. Cc, 89.75.Kd 



I. INTRODUCTION 

Network theory has proven a powerful framework for 
studying the effects of randomness and heterogeneity on 
the dynamics of interacting agents with non-trivial con- 
nectivity patterns [1}. One of the most important ap- 
plications of this work is to the spread of infectious dis- 
eases among human populations, where the interaction 
structure is highly complex, showing salient features such 
as power-law degree distributions, small average path 
lengths, and modularity [2, 3 . Various models have been 
proposed, primarily with random graph configuration, 
that incorporate these complex features while remain- 
ing theoretically tractable. Within the context of disease 
dynamics, graph nodes are generally taken to represent 
individuals, and edges to represent interactions between 
them, through which infection can spread. Both de- 
terministic and stochastic infection dynamics have been 
studied on networks as well as bond percolation for the 
associated branching process O El [TQHT2] . How vari- 
ous thermodynamic quantities of interest - such as the 
steady state incidence, the epidemic (percolation) thresh- 
old, and the distribution of small outbreak sizes - depend 
upon network topology is of great interest. 

Often these approaches disregard the multiscale orga- 
nization of many real systems, in which agents can be 
most naturally thought of as partitioned into densely con- 
nected communities with sparser coupling among neigh- 
boring communities. In some cases, it may be useful 
to conceptualize the topology as a network of networks. 



where agent-to- agent interactions and community-to- 
community interactions are both useful representations 
depending on the scale of resolution [8] . The latter has 
been successfully developed in ecology, with a network of 
interconnected populations referred to as a "metapopula- 
tion" [17l [18] . This framework is very useful in studying 
large scale propagation of diseases where most infection 
transmission occurs in localized regions, but can be trans- 
ported on larger scales by the mobility of individuals, 
traveling among population centers [29 . However, most 
metapopulation models assume that populations are fully 
mixed, with no inherent complexity in the connectivity 
between agents. Much less understood is how the mul- 
tiscale structure of agent interactions affects the larger 
scale propagation of infectious processes through inter- 
connected networks [24^ 26^ . 

In this paper, we expand on a possible avenue for 
addressing this question using a multitype generaliza- 
tion of random graphs with simple, meta-level topology 
O [9] , and construct a dynamical mean-field theory for 
the SIR infection model in multitype configuration model 
networks. Putting these together, we analyze the aver- 
age infection dynamics and propagating front profile on 
a simple metapopulation composed of coupled popula- 
tion centers on a one-dimensional lattice and calculate 
the phenomenological transport properties of the system 
as functions of the underlying network's degree distri- 
butions. Our results are compared to stochastic simu- 
lations of the infection kinetics on various networks and 
found to be in good agreement in the thermodynamic 



limit. Broadly, we present this work as an illustration of 
how well-developed ideas from different areas of statisti- 
cal physics and ecology can be naturally combined. 



By construction, this framework lacks two-point correla- 
tions, in which the excess degree distributions depend on 
the degrees of both nodes sharing an edge [3] . 



II. MULTITYPE CONFIGURATION MODEL 
NETWORKS 



III. VOLZ-MILLER MEAN-FIELD SIR IN 
MULTITYPE NETWORKS 



In order to incorporate relevant node attribute infor- 
mation into our network models, (generically applicable 
for such things as age, sex, ethnicity, and place of res- 
idence), we use a generalization of configuration model 
random graphs, wherein nodes are assigned a type from 
an arbitrary set of M possible types and a degree to 
each type from an arbitrary joint distribution for de- 
gree types, Pi(/ci, /c2, .../cm) = Pi{k)^ with degree kj de- 
noting the number of connections to nodes of type j 
[H [HI 0. Additionally, nodes of type i occupy a frac- 
tion of the total network Wi^ where XI i '^i — 1- ^^^~ 
lowing the configuration model prescription, we consider 
graphs chosen uniformly at random from the ensemble of 
possible graphs with the prescribed degree distributions 
and self-consistent edge constraint: Wi^^kjPi{k) = 



w. 
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From this formalism, a variety of quantities can be de- 
scribed compactly using generating functions p3 17 . The 
generating function for the probability of a randomly se- 
lected node of type i to have degree ^, is given by 



M 



G,(f)=^p,(fc)n^ 



(1) 
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written as a power series in x, an auxiliary variable de- 
fined over the unit interval, with expansion coefficients 
equal to the respective probabilities. Moments of the 
degree distributions can be represented simply as deriva- 
tives of the corresponding generating function. For ex- 
ample, the average degree of a type i node to a type j 
node is 



Y,kjPi{k) = d,^Gi{x)\~^ = {k^), 



(2) 



Since node interactions occur along edges, an important 
quantity in network models is the excess degree: the num- 
ber of neighbors a node has which can be reached by se- 
lecting a randomly chosen edge, and not including the 
neighbor on the end of the selected edge. For a multiype 
configuration model network, the probability that a ran- 
domly chosen edge from a type i node leads to a type j 
node with degree k is proportional to kiPj{k), and thus 
the probability for the corresponding excess degree is gen- 
erated by dxiGj{x)/dxiGj{x)\:^ [9 , with average degree 
to type / nodes. 
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(3) 



In this report we consider simple dynamics for dis- 
ease spread: the Susceptible-Infected- Recovered (SIR) 
model, wherein each individual is assigned a disease state, 
Y G {S^ /, i?}, and may undergo reactions to other states 
depending on its state and the state of its neighbors. In 
this model, if a node of type i is susceptible and has a 
single infected neighbor of type j, then it will change 
its state to infected with a constant probability per unit 
time Pji. Likewise, an infected node of type i will recover 
with a constant probability per unit time 7^. Since the 
underlying dynamics is a continuous time Markov pro- 
cess, a complete analysis would describe the full prob- 
ability distribution for all system trajectories. However 
for our purposes, it will be sufficient to focus on the be- 
havior of extensive outbreaks (i.e., those which scale with 
the system size) , the average dynamics of which, can be 
derived in the limit when the number of nodes tends to 
infinity, by generalizing a mean-field technique for single 
type networks, developed by Volz and Miller, to multi- 
type networks. Below, we follow the basic structure of 
the derivations presented in [13, 14 . 

In the thermodynamic limit, configuration model ran- 
dom graphs are locally tree-like [4] , which by construction 
allows them to satisfy many of the generic criteria for 
the applicability of mean-field theory assumptions [15]. 
In our case, we assume that nodes are differentiated by 
their degree and disease state alone and that suscepti- 
ble nodes feel a uniform force of infection along every 
edge, related to the average number of edges connecting 
susceptible and infected nodes at any given time in the 
network: a Curie- Weiss type approximation [16]. Fur- 
thermore, from the perspective of susceptible nodes, all 
infection attempts along different edges can be treated as 
uncorrelated - a consequence of the local tree-like prop- 
erty [H ini [13]- and thus we assume that the states of 
neighbors of susceptible nodes are effectively indepen- 
dent. 

Let the probability that a node of type j has not trans- 
mitted the infection to a node of type i along a randomly 
chosen i—j edge, be Oij. This quantity is interpretable as 
the complement of the average cumulative hazard func- 
tion along such edges. Given 9ij^ it follows that the frac- 
tion of susceptible nodes of type i at time t is 

M 

s,{t) = ^P^(^) n^'/ W = G^{e,^(t),e,2it),...e,Mit)) 
k i=i 



G,(e,(t)). 



(4) 



The fractions of infected and recovered nodes of type i 
follow from probability conservation. Si -\- li -\- Ri = 1, 



and a constant recovery rate for infected nodes 7^: 



dli 
~dt 
dRi 

~dt~ 



rif) 



dt 



with the total fraction of susceptible nodes 



(5) 



(6) 



The central probability and order parameter, Oij^ can 
be subdivided into three compartments depending on the 
disease state of the terminal node j, 



9^. 



9'- 
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(7) 



and its dynamics determined by tracking the fluxes 
among these compartments. Since 6 can only change 
when an infected node transmits the disease, the rate 
at which Oij changes is equal to the rate at which a 
corresponding neighbor infects, and therefore dOij = 
—pji^ijdt. Similarly, since 6^ can only change if an in- 
fected node recovers, the rate at which 0^ changes is 
equal to the rate at which a corresponding neighbor re- 
covers, and thus d6(j = jjO(jdt. Lastly, the probability 
that a type j neighbor of a type i node has not trans- 
mitted and is susceptible. Of,, is simply the probabil- 
ity that the corresponding neighbor is susceptible. Be- 
cause this neighbor could not have been infected along 
any of its other edges and has excess degree distribu- 
tion generated by dxiGj{x)/dxiGj{x)\:^ , it follows that 
6fj = dxiGj{x)\g-*/dxiGj{x)\:^. Combining the latter 

with the two flux relations and the initial conditions (JTl), 
Oij{0) = 1 and i9^^(0) = 0, we find 



dOij 

~dr 



PJ^ 



'dx,Gj{x)\^^ 
9x,Gj{x)\^ 



+ 7.-(l-%). 



(8) 



These M^, first-order and coupled DDEs, £ = £(£), de- 
fine the full system's approximate mean dynamics, and 
form the basis of our subsequent analysis. For a more de- 
tailed derivation of the analogous results for the special 
case of a single type network, see [13, 14 . 

The steady state is given by the fixed point of (Isl), 



(i-T,,) 






(9) 



which upon substitution into Q, gives the cumulative 
infection, P = 1 — 5, at equilibrium (i.e., the final epi- 
demic size), with Tji = j3ji/ {13 ji + 7^) the corresponding 
bond percolation probability, or transmissibility [9 . This 
can have a non-trivial solution corresponding to the ex- 
istence of extensive outbreaks, if the disease- free state. 



Oij = 1 y{i,j), is unstable. The threshold or phase tran- 
sition, which signifies the region in parameter space that 
separates the epidemic and non-epidemic phases, can be 
obtained through a stability analysis of the disease-free 
state, where the eigenvalue of the Jacobian for Q with 
the largest real part, is real and vanishes when 



det{K - Z) = 0, with 

^^(i,j)(/c,o = TjiSjk {ki)i_j 



(10) 



an M^xM^ matrix [32 . Similar results for the equilib- 
rium properties are derivable from a multitype bond per- 
colation approach [9J. 



IV. FRAMEWORK FOR MULTISCALE 
NETWORKS 

Of interest to us are systems where type structure adds 
an additional scale of relevant topology, and not just de- 
mographic complexity [3 [9]. For instance, we can apply 
the multitype network formalism to a simple model for 
a metapopulation by affiliating population centers with 
node types and coupling among populations with edges 
connecting their constituent nodes. In this way, a com- 
plex topology can be encoded on a micro-scale with a 
macro-scale adjacency matrix, A, describing which pop- 
ulations are directly connected through node interactions 
[8 . We envisage example systems where A describes the 
connectivity among urban centers, such as cities, towns, 
or villages, facilitated by roads or airlines. By concep- 
tualizing the topology in this manner, we can study the 
phenomenology of infection propagation among popula- 
tion centers and describe how the propagation proper- 
ties depend on the underlying connectivity patterns. A 
schematic is shown in Fig{l]-(a) for a simple system with 
the pertinent structure. More broadly, we advance this 
approach as an avenue for combining the frameworks of 
network theory, metapopulations, and front propagation, 
which will be particularly useful if the interaction topol- 
ogy is coherent after some level of coarse graining. 



V. 1-D LATTICE METAPOPULATION 
DYNAMICS 

To illustrate this approach, we consider a special case 
of the above where the macro-scale topology is an infi- 
nite one-dimensional lattice, M ^ 00, in which agents 
interact with other agents of the same type and agents of 
neighboring types, A^j = (^j,n + ^j,n+i + ^j,n-i) [29^. If 
infection is started at a single site (e.g., site 0) in a fully 
susceptible system, a strict directionality applies: in or- 
der for site i to be infected, sites i — 1, i — 2, ... must be 
infected first. In such a case, we expect a well-defined in- 
fection front to propagate through the lattice. In keeping 
with the above, we focus on an effective force of infection 
model among populations with static configuration net- 




{Onn.Onn+l.Onn-l) = (^^(t), ^ W, ^n W) , where (0), 

(+), and (— ) denote internal, right-external, and left- 
external edges, at the corresponding site. For simplicity, 
homogeneity is assumed, with /3, w^ 7, and G all uniform 
- reducing the field equations to 



dr 



{i-et) + T 






(11) 



where the time, r, is measured in units of l/(/3 + 7), 
and the subscript in G denotes a partial derivative with 
respect to the corresponding variable. For edge number 
consistency, G+(l) = G_(l), but in general we allow for 
other asymmetries in the degree distributions. 



A. Dispersion Relation and Transport Coefficients 



(b) 



To understand the spatio-temporal dynamics (11), 




internal edge • susceptible 
. external edge :«-4^, 



FIG. 1: (a) A schematic of SIR dynamics on a metapopu- 
lation, where infection spreads along edges connecting nodes 
of various types at the finest scale (shown with integer, pop- 
ulation labels) , and the macro-scale topology identifies which 
populations are directly connected through agent interactions. 
(b) A particular example of this framework, in which the 
macro-scale topology takes the form of a one-dimensional lat- 
tice. In Sec|V| we focus on a simple case with configura- 
tion model construction, where each site has an identical de- 
gree distribution, specifying the probability of having a given 
number of internal (0), right (+), and left (-) external edges 
(shown above with labels for site i). 



works having prescribed degree distributions - a gener- 
alization of the paradigmatic, spatial SIR model in one 
dimension, where the assumption of well mixed popula- 
tions is relaxed to include complexity in agent interac- 
tions [l8l[T9]. A schematic is shown in Fig{l]-(b). 

Since each node has three edge types, the mean equa- 
tions of motion describe a three-component field, 9n{t) = 



we first quantify how perturbations away from the 
unstable state propagate by linearizing the dynamics 
around the disease-free equilibrium, On{t) = 1 — e'^(t), 
and decoupling the perturbations into basis modes us- 
ing the Inverse Discrete Fourier Transform, en{t) = 
^E!l"^^e;(t)e^(2^^^/^), (IDFT). The dispersion rela- 



M 



tion can be found by substituting the IDFT into (11) 



and using the basis properties of orthogonality and com- 
pleteness. In the limit M ^ 00 the site perturbations 
approach the integral, e^it) = ^ /q ^ e{k)e^^^^~^^^^^^dk. 
With this prescription, we find the dispersion relation 
takes the form of a cubic, characteristic equation 
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which for convenience, is written in terms of s and g, 
where u = is and k = iq. Interestingly, this method 
reveals a generalization of the average excess degree ma- 
trix, -K^(O) - whose elements are found by selecting a ran- 
domly chosen edge in a particular direction, and counting 
the average number of reachable neighbors of a particu- 
lar type - for the interconnected network system, Ke{q), 
which incorporates the relative states of adjacent sites on 
the lattice for each mode q. We expect this operator to 
emerge in similar problems on interconnected networks. 
Combining the above with the behavior of infection 
near the phase transition, where there is no exponential 



growth in time and each site has the same field value: 
5^0 and q ^ ([l2|, we find a simple condition for the 
critical transmissibility Tc'. 



1 



A^(0) ' 



(13) 



where X^{q) is the maximum eigenvalue of Ke{q)^ with 
A^(0) corresponding to -^(0). Because the addition of 
external edges increases the spreading capacity of the dis- 
ease, the critical transmissibility in the coupled system 
is less than the uncoupled case, implying that transport- 
mediated infections from neighboring sites can sustain 
epidemics even when individual populations on their own 
cannot [241126] . 

Also from the dispersion relation, we can find the 
asymptotic transport coefficients for rightward mov- 
ing disturbances by making a standard saddle-point- 
approximation of the perturbations' integral representa- 
tion in fourier space: expanding the integrand around 
its dominant contribution, /c*, in the co- moving frame. 



i{kn-^{k)t) ik^ it[kv*-uj{k*)- 

e ~e e 



du! 
' dk 



(k-k*) 



-it{k-k*f d'^OJ 
2 d/e2 



and taking the infinite time limit while enforcing approxi- 
mate constancy with no exponential growth and ^ finite - 
where v* is the asymptotic speed at which perturbations 
to the unstable state propagate [21^. This procedure un- 
covers an exponential moving pulse for the leading-edge 
of the infection profile with a diffusive correction [21] : 



e-Q^^e-^'/"^^*^ 



(14) 



where g*, 'U*, and I)* satisfy the saddle-point relations: 



ds 
dq 



siq*) 



d\l 



q" 



dq 



-l + T\Uq*) 



and D* = 



2dq^ 



2 dq^ 



(15) 



(16) 



giving a transcendental equation for (7*. When the av- 
erage excess degree matrix is irreducible (the domain of 
interest to us), the dominant growth exponent for each 
q is real and corresponds to a uniquely positive eigen- 
vector (12) [32 , and thus we expect the same selected 



velocity for all fields [23 . Approximately, the fields prop- 
agate in this regime with proportions 1 — e(t)Q((7*, 5*), 



where Q(^*,5*) is the corresponding mode of K^(g*)(12). 
If multiple solutions exist for 1;*, the fastest solution is 
selected [21 . The characteristic wavelength, l/q*, is re- 
lated to the asymptotic size of the front's leading edge, 
and diverges near the phase transition. The diffusion 



coefficient, D*, gives the effective widening of the mean- 
field pulse in the co-moving frame and is proportional to 
the largest finite-size correction to 'U* in the limit where 
the number of nodes at each site tends to infinity. 

In order to uncover the principal dependencies of the 
transport coefficients, we study (fT2|)-(16) near the phase 



transition, where the power series expansion for the dis- 
persion relation is a convenient representation; the latter 
is found by substituting {s{q) + 1)/T = a-^bq-^ |g^ + ... 
into (12), and equating powers in q. When 5* and q* are 



small in the vicinity of Tc, the expansion can be truncated 
at low order, giving a Fisher-Kolmogorov-like dispersion 
relation with the approximate scaling 



s* ^ (rA^(o) - 1) 
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(Al(0)-At(0))(Ai(0)-A3^(0)) 



(17) 



where A2(0) and A^O) are the subdominant eigenvalues 
of Ke{0). In this regime, we find an effective reaction- 
diffusion behavior with the generic dependence of the 
shape and speed of the propagating front's leading edge 
on the reproductive number TA^(O) (a product of the 
spreading capacity along edges and the magnitude of 
topological ffuctuations) and on the normalized diffusion 
coefficient 5: measuring the relative strength of connec- 
tion between lattice sites (17). We see that the effective 



reaction rate is equal to the distance from the phase tran- 
sition, TA^(O) — 1, and that all coefficients grow from zero 
with this distance, except for D* which varies discontin- 
uously through Tc. Furthermore, the normalized diffu- 
sion coefficient increases from zero with (k-ko) (kok-^) 
and (k-k-^) - the correlation moments of the degree dis- 
tribution which encode the propensity for transport from 
the z T 1 site to the i ± 1 site (both of which cannot be 
zero, otherwise epidemics are locally confined), and with 
the viability of subdominant modes to support growth. 
In general, we find that as S increases: v* and D* in- 
crease, g* decreases, and 5* remains constant, implying 
faster transport and greater similarity among sites, as 
more edge- type pairs allow for traversing the lattice, but 
with little change in the growth exponent. 

The above demonstrates the typical trend for these 
models, that the front dynamics is strongly infiuenced 
by the joint degree distribution's second moments (i.e., 
the relevant excess degree properties are generally ampli- 
fied by correlation among degree types and degree hetero- 
geneity). For example, in analogy with the single network 
case, fast transport can be achieved with the presence of 
a few nodes with large internal and external degrees, or 
"transport hubs" , even if the average degrees in the net- 
work are small ^ ,3^ . 



B. Simple Mixing Example 

Additional understanding of the basic form of 
the transport coefficients is gained by looking at 
a special case of the micro-scale degree distribu- 
tion, where the generating function takes the form 

G ( fxQ -\ — Y~ (^+ + ^-) ) 7 with total degree described by 
G, and a given edge connecting nodes of the same site 
with probability /, and nodes of left and right neighbor- 
ing sites with equal probability (1 — /) /2, where 1 — / is 
an effective mixing parameter among populations. With 
this prescription, the critical transmissibility is reduced 
to the inverse of the total-edge excess degree, Tc = QnW\ , 
and the normalized diffusion coefficient, to the fraction of 
external edges in a each direction, S = {1 — f)/2. More- 
over, the dispersion relation takes the instructive form 



siq) = -1 



TG"{1) 



(/ + (l-/)cosh(g)) , (18) 



where 5 + 1 is given by the basic reproductive number 
multiplied by the average relative incidence, e~^^^, on 
the end of a randomly selected edge - illustrating the 
intuitive generalization of the single network case, where 
different edge types are more and less likely to connect 
to infected nodes depending on their place in the lattice, 
and thus to contribute to local growth. 

Likewise, from (15) and (16), we find the speed and 
diffusion coefficient. 






and 






2G'(1) 



(19) 



(20) 



where q^ satisfies (15), and v^ is given by the basic re- 



productive number multiplied by the average product of 
relative position and incidence, —Axe~^^^ , on the end 
of a randomly selected edge. Fig. [2] shows the transport 
coefficients, (18), (19), and ( [2Q| ), as functions of T/Tc 
and /, with partial scaling collapse (fl7|) for the corre- 
sponding class of network configurations. The expected 
reaction-diffusion scaling can be observed near the crit- 
ical point, and far away from the critical region, when 
T > Tc, ^*, f$^, ^U^^ and f^^ tend to limit- 
ing curves which depend only on 1 — / , suggesting the 
intuitive asymptotic proportionality to the reproductive 
number. 



C. Pulled Front Classification 

In order to connect the transport properties of the lin- 
ear equations to the full nonlinear system, we refer to the 
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FIG. 2: The scaled transport coefficients for a one- 
dimensional lattice of configuration model networks with ar- 
bitrary total-degree distribution and inter-populatio n mix ing 
parameter 1 — /, shown as functions of the latter (Sec V B ) : q* 
(a), s* (b), i;* (c), and D"" (d). The colored regions mark the 
range of each coefficient, which are bounded by the critical- 
region scaling (T > Tc), and the limiting behavior (T ^ Tc), 
delineated by dashed and solid curves respectively; the for- 
mer are stra ight lines, signifying agreement with the predicted 
scaling (17). Each panel's arrow indicates the direction of in- 
crease in the distance from the phase transition, T/Tc — 1. 



classification of fronts propagating into unstable states, 
which in our system is the fully susceptible metapop- 
ulation lying ahead of the infection front. In general, 
there are two types of deterministic fronts: pulled and 
pushed, with the former having fronts with asymptotic 
speed equal to the linear spreading speed and the lat- 
ter having fronts with asymptotic speed greater than the 
linear spreading speed [21 . Pushed fronts occur because 
nonlinearities in the equations of motion tend to increase 
the growth of perturbations on the unstable state, re- 
sulting in nontrivial front shape dependence of the speed. 
However, in our system all nonlinearities are proportional 
to probability generating functions, ~ G'{9)/G'{1) ( pTJ ), 
which are monotonically increasing over the unit interval. 
Therefore, all nonlinearities tend to increase 6>, and con- 
sequently dampen the growth of infection - a sufficient 
condition for pulled fronts [22 , and thus we anticipate 
fronts in this model to be pulled; this agrees with the 
intuition that epidemic propagation is governed by its 
behavior in a fully susceptible population. In practice, 
the classification has importance for control strategies in 
systems with similar structure, implying that to mitigate 
the spread of infection among populations, efforts should 
be focused on the leading edge of the front, and not on 
larger outbreaks occurring farther behind. 



D. Relaxation Properties 



In addition to quantifying the transport, the front 
speed can be used to extract information about the dy- 
namics away from the unstable state. As shown above, 
the ^-field settles onto a solution with translational sim- 
ilarity, On±x{t) = On{t T -^)i after an initial transient 
period. Behind the leading edge of the front, the behav- 
ior resembles a relaxation to the stable equilibrium ([9]), 
On{t) ^ 1+ f]{t - ^)^0^ f^^-z%n-vH)^ ^Yieie the spa- 
tial rate, |2:*|, is the dominant eigenvalue of the nonlinear 
eigenvalue equation 
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The latter is the analogue of Ke{q) at the stable state, 
which does not depend on the first two moments of the 
degree distribution directly, but on the generating func- 
tion's properties near the equilibrium ([9]). In general, 
the two characteristic spatial rates for this system are 
not equal, |z*| 7^ g*, and when their difference is large, it 
often signifies a significant separation in the time scales 
of growth, 1/5*, and relaxation l/v*|z*|. The latter pro- 
vides an estimate for the amount of time a site is infec- 
tious, with 1/1 2;* I yielding a related estimate for the width 
of the propagating front (i.e., the typical spatial extent 
of an outbreak at a given time). In particular, when the 
front speed is very fast and the degree distribution's sec- 
ond moments are large with the first moments 0(1), we 
find that |2;*| ^ g*, which suggests broad front profiles. 
In this case, propagation and relaxation can be thought 
of as approximately distinct processes. 



VI. COMPARISON WITH STOCHASTIC 
SIMULATIONS 

The above predictions for the mean-field dynamics on 
the one-dimensional metapopulation were compared to 
stochastic simulations of SIR dynamics on random in- 
stances of multi-scale, metapopulation networks, using 
Gillespie's Direct Method [IHl Ell EH] • The graphs were 
constructed using the multitype configuration model by 
first generating a degree sequence from the desired degree 
distribution and then connecting pairs of edge- "stubs" , 
selected uniformly at random [21 [H [9] . An outbreak was 
started by choosing one node from the first lattice site to 
be infected with all others susceptible. Only outbreaks 
which lead to epidemics with 0{N) cumulative infection 




FIG. 3: (a) A comparison between the average cumulative 
infection profile for stochastic simulations of SIR dynamics 
on a one-dimensional lattice of scale-free (blue) and Poisson 
(green) networks and mean- field predictions ( |11| ). Various 
site sizes are shown with different symbols and color shades - 
varying from light to dark for 10^ to 4 x 10^ respectively. Front 
shapes for increasingly large sizes are found to converge to the 
respective mean- field front. The parameters for the underly- 
ing graphs were chosen to be: iT = 100, a = 2 and p = 0.3 for 
the scal e- free, and C = 2.90157 and p = 0.3 for the Poisson 
(Sec VI). A lattice size of 50 sites was used, which was large 



enough to ensure uniformity with the above graph parameters 
and reaction rates /3 = 7 = 1. The arrow indicates the propa- 
gation direction, (b) Stochastic front realizations conditioned 
on the middle lattice site having cumulative infection equal 
to half the equilibrium value (|9| . Averaging over realizations 
produced profiles like those in (a). 



were considered for comparison with mean-field predic- 
tions. In order to ignore fluctuations in the initial tran- 
sients, time was zeroed after the first 100 reactions. 

We are interested in the average shape of the front that 
connects the fully susceptible unstable state lying ahead 
of the infectious wave and the fully recovered (equilib- 
rium) state lying behind it. The average shape was com- 
puted by taking instantaneous "snapshots" of the profile 
for each stochastic realization, conditioned on the mid- 
dle lattice site having cumulative infection equal to half 
the equilibrium value (|9]), and averaging the cumulative 
infection of the other sites over different realizations. In 
general, the "snapshots" did not occur at the same in- 
stant; however, the shifting of different front realizations 
in time, such that they overlapped at a given point, and 
conditionally averaging over the shape, eliminated some 
of the effects of diffusive wandering. The measured fronts 
were compared to the mean-field profiles by integrating 
the lattice equations ( 11 ). A comparison is shown in Fig. 
[3] for two graphs with scale- free and Poisson degree dis- 
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FIG. 4: (a) Convergence of the average velocities, Vp^, to the 
mean-field predictions, v* , for scale-free (blue) and Poisson 
(green) networks (Fig. Isl as functions of the cumulative num- 
ber of infected nodes at each site, shown with fits to the ex- 
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FIG. 5: The average epidemic fraction at each site for a net- 
work with asymmetric generating function (Sec VI) and vary- 
ing transmissibilities. The critical point at which the epidemic 
vanished agrees with the prediction, Tc = 0.115 (13). Each 
site was occupied by 20,000 nodes on a lattice of 100 sites. 



tributions, with generating functions 
Gs.fX^) =Lia (e~^^^xo{l -p + px+)(l -p^px_) 
and 

Gp{x) = exp(C(xo (1 -p + px+) (1 -p^px-) - 1)) , 

where Lic^ is the polylogarithm function with exponent a 
[2] . The parameters for the degree distributions were cho- 
sen such that each network had the same average degree 
and cloning parameter, p (i.e., given a specified inter- 
nal degree distribution, each of a node's internal edges is 
copied to form an external edge with probability p) , but 
with different inherent levels of heterogeneity. 

We see in Fig. [3] that the epidemic front is broader 
in the scale-free network than in the Poisson. This dif- 
ference comes from the much larger front speed of the 
former, which had average excess degrees an order of 
magnitude larger than the latter, |12]) and (15), and the 



relatively similar relaxation times ( |2ip for the two classes 
of networks (implying that the time scale over which a 
site is infectious in each network is roughly the same). 
In the more homogeneous Poisson networks, the front is 
more narrow and propagates through the lattice on the 
same time scales as the local infection dynamics; whereas 
in the scale-free case, the leading edge of the front propa- 
gates quickly through the lattice, followed by a slower re- 
laxation to the stable equilibrium state behind the front. 
This comparison shows that assumptions of homogene- 
ity can drastically underestimate the speed and extent of 
fronts in systems with heterogenous interactions. 

Additionally, the front speed v was numerically esti- 
mated from the average time (rprog) required for the lead- 
ing edge of the front to move forward by one lattice site 



(where the leading edge was defined as that site where 
the incidence first reached a set 0{1) level) and averag- 
ing over such levels; i.e., 1/v = (rprog) ^ once the initial 
spatial variation had decayed. Fig. [4] shows the con- 
vergence of the measured speed from simulations to the 
mean-field prediction for each graph as a function of the 
steady state, cumulative in fecte d population size at every 
lattice site, N = PN (Sec|lll|, with total size N. The 
lines represent fits to the expected scaling of the largest 
finite-size correction for pulled fronts, 'U* 



D*c 



obtained from a general 1/A^ cutoff in the mean- field 
equations [2QH22] : the coefficients are found to be 0{1) 
of the expected scaling. In general, higher order correc- 
tions in N must be calculated from an analysis of the full, 
stochastic system [22 . The very slow convergence in N 
comes from the transport dependence on the linearized 
equations where infinitesimal infection levels apply and 
sensitivity to stochastic effects is highest. This can be 
seen in the fairly large finite-size corrections to the ve- 
locity, particularly for the scale-free network, leading to 
a more narrow conditionally averaged front relative to 
the mean-field, with fewer sites initiated at a given time 
(Fig. §. 

Finally, the average epidemic profile and transmissi- 
bility threshold, (|9| and (13), were compared to simula- 
tions. Fig. |5] plots those comparisons for a system with 
left-right asymmetric generating function 
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on a lattice of 100 sites, with 20,000 nodes on each site. 
Both the epidemic size for various transmissibilities and 
the threshold were found to be in good agreement with 



mean- field predictions, thougii tiie finite-size effects be- 
came increasingly important as the critical region was ap- 
proached, leading to significantly smaller outbreaks near 
the edges of the lattice. 



VII. CONCLUSION 

In this paper we have generalized a mean-field the- 
ory for infection dynamics on multitype networks, and 
used such networks to model multiscale metapopulations. 
Together, this enabled us to explore how macro-scale 
disease propagation depends on micro-scale interaction 
structure. As a necessary first step in this direction, we 
applied the approach to a simple metapopulation model 
for a chain of coupled populations, and derived the trans- 
port properties for infection, including their scaling with 
the disease transmissibility and the statistical properties 
of the underlying network. We also found a threshold for 
the viability of epidemics, and calculated the relaxation 
properties of the propagating front. These were com- 
pared for different network models, with heterogeneous 
networks having considerably higher speeds and broader 
fronts than their homogeneous counterparts - illustrating 
the importance of including complexity in the fine-scale 
topology in order to accurately capture transport phe- 



nomenology. 

Various extensions of the work presented - both in 
terms of analyses carried out and systems studied - could 
be considered. We have addressed here only the aver- 
age dynamics of the one-dimensional, homogeneous sys- 
tem, without any description of finite-size fluctuations, 
or consideration of the dynamics in higher-dimensional 
generalizations. Greater complexity could be introduced 
through the spatio-temporal dependence of network pa- 
rameters, and/or more general network configurations 
[25] . An interesting extension of the model discussed here 
would include dynamic contacts between nodes and ex- 
plicit mobility, instead of the assumed time scale separa- 
tion between topology and the overlying process [29ti3T] . 
However, the basic formalism presented here can enable 
one to study such factors and build more realistic models 
for infectious processes in multiscale problems. 
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